其他
宏基因组02. HUMAnN2 --宏基因组代谢通路分析
本系列课程前情回顾
HUMAnN2 简介
HUMAnN2 全称(HMP Unified Metabolic Analysis Network)
是Huttenhower实验室开发的一款从宏基因组测序数据(通常是数百万短DNA/RNA读长)中高效、准确地描述群落中存在和缺失的微生物以及他们的丰度的流程,这一过程被称为功能特征分析,目的是描述微生物群落及其成员的代谢潜力。更广泛地说,功能分析回答了一个问题:“所感兴趣的微生物群落是做什么的(或者有能力做什么)?”。
HUMAnN2 安装
conda install -n bowtie2 humman2 python=2.7
HUMAnN2 的主要工作流示意图
HUMAnN2数据库准备与配置
humann2_databases # 显示可用数据库
wd=/db/humann2
#mkdir -p $wd # 建立下载目录
# 微生物物种核心基因下载 5.37G
humann2_databases --download chocophlan full $wd
# 功能基因diamond索引下载 10.3G
#humann2_databases --download uniref uniref90_diamond $wd
# 设置数据库位置
# 显示参数
humann2_config --print
# 如修改线程数
humann2_config --update run_modes threads 8
humann2_config --update database_folders protein $wd/uniref
humann2_config --update database_folders nucleotide $wd/chocophlan
#将metaphlan2需要的参考数据库链接到humman2安装的目录下
mkdir -p ~/anaconda3/envs/bowtie2/bin/db_v20
ln -s /db/metaphlan2/* ~/anaconda3/envs/bowtie2/bin/db_v20
mkdir -p ~/anaconda3/envs/bowtie2/bin/databases
ln -s /db/metaphlan2/* ~/anaconda3/envs/bowtie2/bin/databases
#激活kneaddata安装的环境
source activate bowtie2
#单独运行示例
humann2 --input ~/meta/raw_data.dir/concat/10_19_.fq --output temp --threads 8
#并行运行示例
parallel -j 12 'humann2 \
--input {} --output temp' \
::: ~/meta/raw_data.dir/concat/*.fq
物种组成分析
mkdir -p meta/3metagenome/metaphlan2 && cd 3metagenome
# 样品结果合并
merge_metaphlan_tables.py temp/*_humann2_temp/*_metaphlan_bugs_list.tsv | sed 's/_metaphlan_bugs_list//g' > metaphlan2/taxonomy.tsv
# 转换为spf格式方便stamp分析
metaphlan_to_stamp.pl metaphlan2/taxonomy.tsv > metaphlan2/taxonomy.spf
# 5. 物种组成分析和可视化进阶
# 绘制热图
metaphlan_hclust_heatmap.py --in metaphlan2/taxonomy.tsv \
--out metaphlan2/heatmap.pdf -c bbcry --top 25 --minv 0.1 -s log
# c设置颜色方案,top设置物种数量,minv最小相对丰度,s标准化方法,log为取10为底对数,文件名结尾可先pdf/png/svg三种图片格式。更多说明详见 metaphlan_hclust_heatmap.py -h
# GraPhlAn图
# metaphlan2 to graphlan
export2graphlan.py --skip_rows 1,2 -i metaphlan2/taxonomy.tsv\
--tree temp.bak/merged_abundance.tree.txt \
--annotation temp.bak/merged_abundance.annot.txt \
--most_abundant 100 --abundance_threshold 1 --least_biomarkers 10 \
--annotations 5,6 --external_annotations 7 --min_clade_size 1
##--skip_rows 输入文件中略过的行数
##-i 输入文件
##--tree 输出的树文件用于GraPhlAn
##--annotation 输出的注释文件文件
##--most_abundant 当只提供Lefse_input时,指定突出显示多少个clades。由于生物标记物缺失,它们将从最丰富的生物标记中选出。缺省值为10。
##--abundance_threshold 为要注释的类设置最小丰度值。缺省值为20.0
# graphlan annotation
graphlan_annotate.py --annot temp.bak/merged_abundance.annot.txt temp.bak/merged_abundance.tree.txt temp.bak/merged_abundance.xml
##--annot注释文件
# output PDF figure, annoat and legend
graphlan.py temp.bak/merged_abundance.xml ./metaphlan2/graphlan.pdf --external_legends --dpi 300
本次学习到此为止,谢谢关注,不足之处请多指正。